Let us set some global options for all code chunks in this
document.
Utilitary functions
# For each m and beta, this function returns c_m/b_{m+1} and the roots of rb and rc
my.get.roots <- function(order, beta) {
mt <- get(paste0("m", order, "table"))
rb <- rep(0, order + 1)
rc <- rep(0, order)
if(order == 1) {
rc = approx(mt$beta, mt[[paste0("rc")]], beta)$y
} else {
rc = sapply(1:order, function(i) {
approx(mt$beta, mt[[paste0("rc.", i)]], beta)$y
})
}
rb = sapply(1:(order+1), function(i) {
approx(mt$beta, mt[[paste0("rb.", i)]], xout = beta)$y
})
factor = approx(mt$beta, mt$factor, xout = beta)$y
return(list(rb = rb, rc = rc, factor = factor))
}
# Given the roots, return the polynomial coefficients in increasing order like a+bx+cx^2+...
# So if the output is c(6, -5, 1), that means 6 - 5x + x^2
poly_from_roots <- function(roots) {
coef <- 1
for (r in roots) {
coef <- convolve(coef, c(1, -r), type = "open")
}
return(coef) # returns in increasing order like a+bx+cx^2+...
}
# Given a set of complex roots, this function groups complex roots into conjugate pairs
get_conjugate_pairs <- function(roots) {
used <- rep(FALSE, length(roots))
pairs <- list()
for (i in seq_along(roots)) {
if (!used[i]) {
conj_root <- Conj(roots[i])
# Find its conjugate (within tolerance)
j <- which(!used & abs(Re(roots) - Re(conj_root)) < 1e-8 & abs(Im(roots) - Im(conj_root)) < 1e-8)
j <- j[j != i]
if (length(j) > 0) {
used[c(i, j[1])] <- TRUE
pairs[[length(pairs) + 1]] <- c(roots[i], roots[j[1]])
}
}
}
return(pairs)
}
# This function builds the polynomial c_m\b_{m+1} prod_{i=1}^m (1-r_{1i}x) + tau \prod_{j=1}^{m+1} (1-r_{2j}x)
# For example compute_sum_poly(0, rep(0,2), c(2,-3,1),1) = c(6,-7,0,1) because (1-2x)(1+3x)(1-x) = 6x^3 - 7x^2 + 0x + 1
compute_sum_poly <- function(factor, pr_roots, pl_roots, cte) {
pr_coef <- c(0, poly_from_roots(pr_roots)) # in decreasing order like x^n+bx^(n-1)+cx^(n-2)+...
pl_coef <- poly_from_roots(pl_roots) # in decreasing order like x^n+bx^(n-1)+cx^(n-2)+...
pr_plus_pl_coef <- factor * pr_coef + cte * pl_coef
return(pr_plus_pl_coef) # in decreasing order like x^n+bx^(n-1)+cx^(n-2)+...
}
compute_real_roots_and_complex_coef <- function(pr_plus_pl_coef) {
pr_plus_pl_roots <- polyroot(rev(pr_plus_pl_coef)) # polyroot expects the coefficients in increasing order
real_roots <- Re(pr_plus_pl_roots[abs(Im(pr_plus_pl_roots)) < 1e-8])
complex_roots <- pr_plus_pl_roots[abs(Im(pr_plus_pl_roots)) >= 1e-8]
complex_poly_coefs <- list()
if (length(complex_roots) > 0) {
pairs <- get_conjugate_pairs(complex_roots)
for (pair in pairs) {
coef <- rev(Re(poly_from_roots(c(pair[1], pair[2])))) # returns in x^2 + bx + c order
complex_poly_coefs[[length(complex_poly_coefs) + 1]] <- round(coef, 12)
}
} else {
complex_poly_coefs <- list() # No complex roots
}
return(list(new_factor = pr_plus_pl_coef[1],
pr_plus_pl_roots = pr_plus_pl_roots,
real_roots = real_roots,
complex_poly_coefs = complex_poly_coefs))
}
my.fractional.operators <- function(L, # kappa^2C + G
beta,
C,
scale.factor, # kappa^2
m = 1,
time_step) {
C <- Matrix::Diagonal(dim(C)[1], rowSums(C)) # lumped
Ci <- Matrix::Diagonal(dim(C)[1], 1 / rowSums(C)) # lumped
I <- Matrix::Diagonal(dim(C)[1])
L <- L / scale.factor # C + G/kappa^2
LCi <- L %*% Ci
roots <- my.get.roots(m, beta)
Pl.roots <- roots$rb
Pr.roots <- roots$rc
factor <- roots$factor
new_factor_and_roots <- compute_real_roots_and_complex_coef(compute_sum_poly(factor, Pr.roots, Pl.roots, time_step))
new_real_roots <- new_factor_and_roots$real_roots
new_factor <- new_factor_and_roots$new_factor
complex_poly_coefs <- new_factor_and_roots$complex_poly_coefs
# Fill Pr_plus_Pl.factors
Pr_plus_Pl.factors <- list()
if (length(new_real_roots) >= 1) {
for (i in 1:length(new_real_roots)) {
Pr_plus_Pl.factors[[i]] <- LCi - new_real_roots[i] * I
}
}
print(length(Pr_plus_Pl.factors))
length_now <- length(Pr_plus_Pl.factors)
if(length(complex_poly_coefs) >= 1) {
for (i in 1:length(complex_poly_coefs)) {
}
Pr_plus_Pl.factors[[length_now + i]] <- LCi %*% LCi + complex_poly_coefs[[i]][2] * LCi + complex_poly_coefs[[i]][3] * I
}
# Fill Pr.factors
Pr.factors <- list()
Pr.factors[[1]] <- I - LCi * roots$rc[1]
if (length(roots$rc) > 1) {
for (i in 2:length(roots$rc)) {
Pr.factors[[i]] <- I - LCi * roots$rc[i]
}
}
return(list(
Ci = Ci, C = C, LCi = LCi, L = L, m = m, beta = beta,
factor = factor,
new_factor = new_factor,
Pr.factors = Pr.factors,
Pr_plus_Pl.factors = Pr_plus_Pl.factors
))
}
my.solver <- function(obj, v){
Pr.factors <- obj$Pr.factors
Pr_plus_Pl.factors <- obj$Pr_plus_Pl.factors
m <- obj$m
C <- obj$C
Ci <- obj$Ci
factor <- obj$factor
new_factor <- obj$new_factor
if (m==1){
temp <- solve(Pr_plus_Pl.factors[[2]],
Pr.factors[[1]] %*% solve(
Pr_plus_Pl.factors[[1]],
C%*%v))
} else if (m==2){
temp <- solve(Pr_plus_Pl.factors[[3]],
Pr.factors[[2]] %*% solve(
Pr_plus_Pl.factors[[2]],
Pr.factors[[1]] %*% solve(
Pr_plus_Pl.factors[[1]],
C%*%v)))
} else if (m==3){
temp <- solve(Pr_plus_Pl.factors[[4]],
Pr.factors[[3]] %*% solve(
Pr_plus_Pl.factors[[3]],
Pr.factors[[2]] %*% solve(
Pr_plus_Pl.factors[[2]],
Pr.factors[[1]] %*% solve(
Pr_plus_Pl.factors[[1]],
C%*%v))))
} else if (m==4){
temp <- solve(Pr_plus_Pl.factors[[5]],
Pr.factors[[4]] %*% solve(
Pr_plus_Pl.factors[[4]],
Pr.factors[[3]] %*% solve(
Pr_plus_Pl.factors[[3]],
Pr.factors[[2]] %*% solve(
Pr_plus_Pl.factors[[2]],
Pr.factors[[1]] %*% solve(
Pr_plus_Pl.factors[[1]],
C%*%v)))))
}
return((factor/new_factor) * Ci %*% temp)
}
We want to solve the fractional diffusion equation \[\begin{equation}
\label{eq:maineq}
\partial_t u+(\kappa^2-\Delta_\Gamma)^{\alpha/2} u=f \text { on }
\Gamma \times(0, T), \quad u(0)=u_0 \text { on } \Gamma,
\end{equation}\] where \(u\)
satisfies the Kirchhoff vertex conditions \[\begin{equation}
\label{eq:Kcond}
\left\{\phi\in C(\Gamma)\;\Big|\; \forall v\in V:
\sum_{e\in\mathcal{E}_v}\partial_e \phi(v)=0 \right\}
\end{equation}\]
If \(f=0\), then the solution is
given by \[\begin{equation}
\label{eq:sol_reprentation}
u(s,t) =
\displaystyle\sum_{j\in\mathbb{N}}e^{-\lambda^{\alpha/2}_jt}\left(u_0,
e_j\right)_{L_2(\Gamma)}e_j(s).
\end{equation}\]
# Function to build a tadpole graph and create a mesh
gets_graph_tadpole <- function(h){
edge1 <- rbind(c(0,0),c(1,0))
theta <- seq(from=-pi,to=pi,length.out = 100)
edge2 <- cbind(1+1/pi+cos(theta)/pi,sin(theta)/pi)
edges = list(edge1, edge2)
graph <- metric_graph$new(edges = edges)
graph$build_mesh(h = h)
return(graph)
}
Let \(\Gamma_T =
(\mathcal{V},\mathcal{E})\) characterize the tadpole graph with
\(\mathcal{V}= \{v_1,v_2\}\) and \(\mathcal{E}= \{e_1,e_2\}\). The left edge
\(e_1\) has length 1 and the circular
edge \(e_2\) has length 2. As discussed
before, a point on \(e_1\) is
parameterized via \(s=\left(e_1,
t\right)\) for \(t \in[0,1]\)
and a point on \(e_2\) via \(s=\left(e_2, t\right)\) for \(t\in[0,2]\). One can verify that \(-\Delta_\Gamma\) has eigenvalues \(0,\left\{(i \pi / 2)^2\right\}_{i \in
\mathbb{N}}\) and \(\left\{(i \pi /
2)^2\right\}_{2 i \in \mathbb{N}}\) with corresponding
eigenfunctions \(\phi_0\), \(\left\{\phi_i\right\}_{i \in \mathbb{N}}\),
and \(\left\{\psi_i\right\}_{2 i \in
\mathbb{N}}\) given by \(\phi_0(s)=1 /
\sqrt{3}\) and \[\begin{equation*}
\phi_i(s)=C_{\phi, i}\begin{cases}
-2 \sin (\frac{i\pi}{2}) \cos (\frac{i \pi t}{2}), & s \in
e_1, \\
\sin (i \pi t / 2), & s \in e_2,
\end{cases},
\quad
\psi_i(s)=\frac{\sqrt{3}}{\sqrt{2}} \begin{cases}
(-1)^{i / 2} \cos (\frac{i \pi t}{2}), & s \in e_1, \\
\cos (\frac{i \pi t}{2}), & s \in e_2,
\end{cases},
\end{equation*}\] where \(C_{\phi,
i}=1\) if \(i\) is even and
\(C_{\phi, i}=1 / \sqrt{3}\) otherwise.
Moreover, these functions form an orthonormal basis for \(L_2(\Gamma_T)\).
# Function to compute the eigenfunctions
tadpole.eig <- function(k,graph){
x1 <- c(0,graph$get_edge_lengths()[1]*graph$mesh$PtE[graph$mesh$PtE[,1]==1,2])
x2 <- c(0,graph$get_edge_lengths()[2]*graph$mesh$PtE[graph$mesh$PtE[,1]==2,2])
if(k==0){
f.e1 <- rep(1,length(x1))
f.e2 <- rep(1,length(x2))
f1 = c(f.e1[1],f.e2[1],f.e1[-1], f.e2[-1])
f = list(phi=f1/sqrt(3))
} else {
f.e1 <- -2*sin(pi*k*1/2)*cos(pi*k*x1/2)
f.e2 <- sin(pi*k*x2/2)
f1 = c(f.e1[1],f.e2[1],f.e1[-1], f.e2[-1])
if((k %% 2)==1){
f = list(phi=f1/sqrt(3))
} else {
f.e1 <- (-1)^{k/2}*cos(pi*k*x1/2)
f.e2 <- cos(pi*k*x2/2)
f2 = c(f.e1[1],f.e2[1],f.e1[-1],f.e2[-1])
f <- list(phi=f1,psi=f2/sqrt(3/2))
}
}
return(f)
}
Let \(\alpha\in(0,2]\) and \(U_h^\tau\) denote the sequence of
approximations of the solution to the weak form of problem \(\eqref{eq:maineq}\) at each time step on a
mesh indexed by \(h\). Let \(z=0\) and \(U^0_h
= P_hu_0\). For \(k=0,\dots,
N-1\), \(U_h^{k+1}\in V_h\)
solves the following scheme \[\begin{align}
\label{system:fully_discrete_scheme}
\langle\delta U_h^{k+1},\phi\rangle +
\mathfrak{a}(U_h^{k+1},\phi) = \langle f^{k+1},\phi\rangle
,\quad\forall\phi\in V_h,
\end{align}\] where \(f^{k+1} =
\displaystyle\dfrac{1}{\tau}\int_{t_k}^{t^{k+1}}f(t)dt\). The
solution can be represented as \[\begin{align*}
U_h^k(s) = \sum_{j=1}^{N_h}u_j^k\psi_j(s)
\end{align*}\] Replacing this into \(\eqref{system:fully_discrete_scheme}\)
yields the following linear system \[\begin{align*}
\sum_{j=1}^{N_h}u_j^{k+1}[(\psi_j,\psi_i)_{L_2(\Gamma)}+
\tau\mathfrak{a}(\psi_j,\psi_i)] =
\sum_{j=1}^{N_h}u_j^{k}(\psi_j,\psi_i)_{L_2(\Gamma)}+\tau(
f^{k+1},\psi_i)_{L_2(\Gamma)}
\end{align*}\] for \(i = 1,\dots,
N_h\). In matrix notation, \[\begin{align}
\label{diff_eq_discrete}
(C+\tau L^{\alpha/2})U^{k+1} = CU^k+\tau F^{k+1},
\end{align}\] where \(C\) has
entries \(C_{ij} =
(\psi_j,\psi_i)_{L_2(\Gamma)}\), \(L^{\alpha/2}\) has entries \(\mathfrak{a}(\psi_j,\psi_i)\), \(U^k\) has entries \(u_j^k\), and \(F^k\) has entries \(( f^{k},\psi_i)_{L_2(\Gamma)}\). For
simplicity, let \(f=0\) and so we
arrive at
\[\begin{equation}
(C+\tau L^{\alpha/2})U^{k+1} = CU^{k}
\end{equation}\]
Apply \(L^{-\alpha/2}\) to both
sides to get
\[\begin{equation}
(L^{-\alpha/2}C+\tau I)U^{k+1} = L^{-\alpha/2}CU^{k}
\end{equation}\]
Approximate \(L^{-\alpha/2}\) by
\(P_rP_l^{-1}\) to arrive at
\[\begin{equation}
(P_rP_l^{-1}C+\tau I)U^{k+1} = P_rP_l^{-1}CU^{k}
\end{equation}\]
Actually, \(L^{-\alpha/2}\) should
be approximated by \(P_l^{-T}P_r^T\)
(this is like this because we want to use \(P_r\) and \(P_l\) the way they where constructed) and
so we can write
\[\begin{equation}
(P_l^{-T}P_r^TC+\tau I)U^{k+1} = P_l^{-T}P_r^TCU^{k}
\end{equation}\]
From this, we arrive at the scheme
\[\begin{equation}
(P_r^TC+\tau P_\ell^T)U^{k+1} = P_r^TCU^k
\label{eq:scheme}
\end{equation}\] where (here we are using \(P_r\) and \(P_l\) the way they were constructed) \[\begin{equation}
P_r = c_m\prod_{i=1}^m (I-r_{1i}C^{-1}L)\quad\text{and}\quad P_\ell =
b_{m+1}C\prod_{j=1}^{m+1} (I-r_{2j}C^{-1}L)
\end{equation}\] Observe that \[\begin{equation}
P_r^T = c_m\prod_{i=1}^m (I-r_{1i}LC^{-1})\quad\text{and}\quad P_\ell^T
= b_{m+1}\prod_{j=1}^{m+1} (I-r_{2j}LC^{-1})\cdot C
\end{equation}\] since \(L\) and
\(C^{-1}\) are symmetric and the
factors in the product commute (because the factors in the product that
make up a polynomial commute, for example, \((1-ax)(1-bx) = (1-bx)(1-ax)\)). Replacing
these two in our scheme we get \[\begin{equation}
\left(\dfrac{c_m}{b_{m+1}}\prod_{i=1}^m (I-r_{1i}LC^{-1})+\tau
\prod_{j=1}^{m+1} (I-r_{2j}LC^{-1})\right)CU^{k+1} =
\dfrac{c_m}{b_{m+1}}\prod_{i=1}^m (I-r_{1i}LC^{-1})\cdot CU^k
\end{equation}\]
We now need to compute the roots of the polynomial \[\begin{equation}
R(x) = \dfrac{c_m}{b_{m+1}} \prod_{i=1}^m (1-r_{1i}x)+\tau
\prod_{j=1}^{m+1} (1-r_{2j}x)
\end{equation}\] \(R\) has
leading coefficient \(\gamma=\tau (-1)^{m+1}
\prod_{j=1}^{m+1}r_{2j}\) and \(m+1\) roots \(R_k\) for \(k=1,\ldots,m+1\). That is, \[\begin{equation}
R(x) = \tau(-1)^{m+1} \prod_{j=1}^{m+1}r_{2j}\prod_{k=1}^{m+1} (x-R_k)=
\gamma\prod_{k=1}^{m+1} (x-R_k)
\end{equation}\] We can then write our scheme as \[\begin{equation}
\left(\gamma\prod_{k=1}^{m+1} (LC^{-1}-R_kI)\right)CU^{k+1} =
\dfrac{c_m}{b_{m+1}}\prod_{i=1}^m (I-r_{1i}LC^{-1})\cdot CU^k
\end{equation}\] That is, \[\begin{equation}
U^{k+1} = \dfrac{c_m}{b_{m+1}}\dfrac{1}{\gamma}C^{-1}\prod_{k=1}^{m+1}
(LC^{-1}-R_kI)^{-1}\prod_{i=1}^m (I-r_{1i}LC^{-1})\cdot CU^k
\end{equation}\]
Roots \(r_{1i}\) and \(r_{2j}\) are real but \(R_k\) can be complex. In that case, we
consider the quadratic terms. For example, if \(m=3\), and \(R\) has one real root \(r\) and two complex conjugate roots \(z\) and \(\bar{z}\) such that \((x-z)(x-\bar{z}) = x^2+ax+b\), we will work
with \[\begin{equation}
R(x) = \gamma(x^2+ax+b)(x-r)
\end{equation}\] instead of \[\begin{equation}
R(x) = \gamma(x-z)(x-\bar{z})(x-r)
\end{equation}\]
h <- 0.01
graph <- gets_graph_tadpole(h = h)
## Starting graph creation...
## LongLat is set to FALSE
## Creating edges...
## Setting edge weights...
## Computing bounding box...
## Setting up edges
## Merging close vertices
## Total construction time: 0.16 secs
## Creating and updating vertices...
## Storing the initial graph...
## Computing the relative positions of the edges...
T_final <- 2
time_step <- 0.01
time_seq <- seq(0, T_final, by = time_step)
# Compute the FEM matrices
graph$compute_fem()
G <- graph$mesh$G
C <- graph$mesh$C
I <- Matrix::Diagonal(nrow(C))
x <- graph$mesh$V[, 1]
y <- graph$mesh$V[, 2]
edge_number <- graph$mesh$VtE[, 1]
pos <- sum(edge_number == 1)+1
order_to_plot <- function(v)return(c(v[1], v[3:pos], v[2], v[(pos+1):length(v)], v[2]))
weights <- graph$mesh$weights
kappa <- 1
alpha <- 1.99 # from 0.5 to 2
m = 4
beta <- alpha/2
L <- kappa^2*C + G
True solution
# Initial condition
U_0 <- 10*exp(-((x-1)^2 + (y)^2))
U_true <- matrix(NA, nrow = nrow(C), ncol = length(time_seq))
U_true[, 1] <- U_0
n_finite <- 100
for (k in 1:(length(time_seq) - 1)) {
aux_k <- rep(0, nrow(C))
for (j in 0:n_finite) {
decay_j <- exp(-time_seq[k+1]*(kappa^2 + (j*pi/2)^2)^(alpha/2))
e_j <- tadpole.eig(j,graph)$phi
aux_k <- aux_k + decay_j*sum(U_0*e_j*weights)*e_j
if (j>0 && (j %% 2 == 0)) {
e_j <- tadpole.eig(j,graph)$psi
aux_k <- aux_k + decay_j*sum(U_0*e_j*weights)*e_j
}
}
U_true[, k + 1] <- aux_k
}
my_op <- my.fractional.operators(L, beta, C, scale.factor = kappa^2, m = m, time_step)
## [1] 5
U_approx1 <- matrix(NA, nrow = nrow(C), ncol = length(time_seq))
U_approx1[, 1] <- U_0
# Time-stepping loop
for (k in 1:(length(time_seq) - 1)) {
U_approx1[, k + 1] <- as.matrix(my.solver(my_op, U_approx1[, k]))
}
{r}
op <- fractional.operators(L, beta, C, scale.factor = kappa^2, m = m)
Pl <- op$Pl
Pr <- op$Pr
Ci <- op$Ci
C <- op$C
{r}
LHS <- t(Pr)%*% C + time_step * t(Pl)
# Initialize U matrix to store solution at each time step
U_approx2 <- matrix(NA, nrow = nrow(C), ncol = length(time_seq))
U_approx2[, 1] <- U_0
# Time-stepping loop
for (k in 1:(length(time_seq) - 1)) {
# Compute the right-hand side for the second equation
RHS <- t(Pr)%*% C %*% U_approx2[, k]
U_approx2[, k + 1] <- as.matrix(solve(LHS, RHS))
}
{r}
Pr.apply.mult <- function(v){return(Pr.mult(op, v))}
Pl.apply.solve <- function(v){return(Pl.solve(op, v))}
PliC <- apply(C, 2, Pl.apply.solve) # PlC^-1
# Precompute the LHS1 matrix
aux <- apply(PliC, 2, Pr.apply.mult) #PrPl^-1C
LHS <- aux + time_step * Matrix::Diagonal(nrow(C))
# Initialize U matrix to store solution at each time step
U_approx2 <- matrix(NA, nrow = nrow(C), ncol = length(time_seq))
U_approx2[, 1] <- U_0
# Time-stepping loop
for (k in 1:(length(time_seq) - 1)) {
# Compute the right-hand side for the second equation
RHS <- aux %*% U_approx2[, k]
U_approx2[, k + 1] <- as.matrix(solve(LHS, RHS))
}
Plot
x <- order_to_plot(x)
y <- order_to_plot(y)
U_approx2 <- U_approx1
max_error_at_each_time1 <- apply((U_true - U_approx1)^2, 2, mean)
max_error_at_each_time2 <- apply((U_true - U_approx2)^2, 2, mean)
max_error_between_both_approx <- apply((U_approx1 - U_approx2)^2, 2, mean)
U_true <- apply(U_true, 2, order_to_plot)
U_approx1 <- apply(U_approx1, 2, order_to_plot)
U_approx2 <- apply(U_approx2, 2, order_to_plot)
# Create interactive plot
fig <- plot_ly()
# Add first line (max_error_at_each_time1)
fig <- fig %>% add_trace(
x = ~time_seq, y = ~max_error_at_each_time1, type = 'scatter', mode = 'lines+markers',
line = list(color = 'red', width = 2),
marker = list(color = 'red', size = 4),
name = "Error True and Approx 1"
)
# Add second line (max_error_at_each_time2)
fig <- fig %>% add_trace(
x = ~time_seq, y = ~max_error_at_each_time2, type = 'scatter', mode = 'lines+markers',
line = list(color = 'blue', width = 2),
marker = list(color = 'blue', size = 4),
name = "Error True and Approx 2"
)
# Add third line (max_error_between_both_approx)
fig <- fig %>% add_trace(
x = ~time_seq, y = ~max_error_between_both_approx, type = 'scatter', mode = 'lines+markers',
line = list(color = 'orange', width = 2),
marker = list(color = 'orange', size = 4),
name = "Error Between Approximations"
)
# Layout
fig <- fig %>% layout(
title = "Error at Each Time Step",
xaxis = list(title = "Time"),
yaxis = list(title = "Error"),
legend = list(x = 0.1, y = 0.9)
)
plot_data <- data.frame(
x = rep(x, times = ncol(U_true)),
y = rep(y, times = ncol(U_true)),
z_true = as.vector(U_true),
z_approx1 = as.vector(U_approx1),
z_approx2 = as.vector(U_approx2),
frame = rep(time_seq, each = length(x))
)
# Compute axis limits
x_range <- range(x)
y_range <- range(y)
z_range <- range(c(U_true, U_approx1, U_approx2))
# Initial plot setup (first frame only)
p <- plot_ly(plot_data, frame = ~frame) %>%
add_trace(
x = ~x, y = ~y, z = ~z_true,
type = "scatter3d", mode = "lines",
name = "True",
line = list(color = "green", width = 2)
) %>%
add_trace(
x = ~x, y = ~y, z = ~z_approx1,
type = "scatter3d", mode = "lines",
name = "Approx 1",
line = list(color = "red", width = 2)
) %>%
add_trace(
x = ~x, y = ~y, z = ~z_approx2,
type = "scatter3d", mode = "lines",
name = "Approx 2",
line = list(color = "blue", width = 2)
) %>%
layout(
scene = list(
xaxis = list(title = "x", range = x_range),
yaxis = list(title = "y", range = y_range),
zaxis = list(title = "Value", range = z_range),
aspectratio = list(x = 2.4, y = 1.2, z = 1.2),
camera = list(
eye = list(x = 1.5, y = 1.5, z = 1), # Adjust the viewpoint
center = list(x = 0, y = 0, z = 0))
),
updatemenus = list(
list(
type = "buttons", showactive = FALSE,
buttons = list(
list(label = "Play", method = "animate",
args = list(NULL, list(frame = list(duration = 100, redraw = TRUE), fromcurrent = TRUE))),
list(label = "Pause", method = "animate",
args = list(NULL, list(mode = "immediate", frame = list(duration = 0), redraw = FALSE)))
)
)
),
title = "Time: 0"
)
# Convert to plotly object with frame info
pb <- plotly_build(p)
# Inject custom titles into each frame
for (i in seq_along(pb$x$frames)) {
t <- time_seq[i]
err <- signif(max_error_between_both_approx[i], 4)
pb$x$frames[[i]]$layout <- list(title = paste0("Time: ", t, " | Error: ", err))
}
LS0tCnRpdGxlOiAiU29sdmluZyBhIHBhcmFib2xpYyBlcXVhdGlvbiIKZGF0ZTogIkNyZWF0ZWQ6IDIwLTA0LTIwMjUuIExhc3QgbW9kaWZpZWQ6IGByIGZvcm1hdChTeXMudGltZSgpLCAnJWQtJW0tJVkuJylgIgpvdXRwdXQ6CiAgaHRtbF9kb2N1bWVudDoKICAgIG1hdGhqYXg6ICJodHRwczovL2Nkbi5qc2RlbGl2ci5uZXQvbnBtL21hdGhqYXhAMy9lczUvdGV4LW1tbC1jaHRtbC5qcyIKICAgIGhpZ2hsaWdodDogcHlnbWVudHMKICAgIHRoZW1lOiBmbGF0bHkKICAgIGNvZGVfZm9sZGluZzogc2hvdyAjIGNsYXNzLnNvdXJjZSA9ICJmb2xkLWhpZGUiIHRvIGhpZGUgY29kZSBhbmQgYWRkIGEgYnV0dG9uIHRvIHNob3cgaXQKICAgIGRmX3ByaW50OiBwYWdlZAogICAgIyB0b2M6IHRydWUKICAgICMgdG9jX2Zsb2F0OgogICAgIyAgIGNvbGxhcHNlZDogdHJ1ZQogICAgIyAgIHNtb290aF9zY3JvbGw6IHRydWUKICAgIG51bWJlcl9zZWN0aW9uczogZmFsc2UKICAgIGZpZ19jYXB0aW9uOiB0cnVlCiAgICBjb2RlX2Rvd25sb2FkOiB0cnVlCmFsd2F5c19hbGxvd19odG1sOiB0cnVlCmJpYmxpb2dyYXBoeTogCiAgLSByZWZlcmVuY2VzLmJpYgogIC0gZ3JhdGVmdWwtcmVmcy5iaWIKaGVhZGVyLWluY2x1ZGVzOgogIC0gXG5ld2NvbW1hbmR7XGFyfXtcbWF0aGJie1J9fQogIC0gXG5ld2NvbW1hbmR7XGxsYXZ9WzFde1xsZWZ0XHsjMVxyaWdodFx9fQogIC0gXG5ld2NvbW1hbmR7XHBhcmV9WzFde1xsZWZ0KCMxXHJpZ2h0KX0KICAtIFxuZXdjb21tYW5ke1xOY2FsfXtcbWF0aGNhbHtOfX0KICAtIFxuZXdjb21tYW5ke1xWY2FsfXtcbWF0aGNhbHtWfX0KICAtIFxuZXdjb21tYW5ke1xFY2FsfXtcbWF0aGNhbHtFfX0KICAtIFxuZXdjb21tYW5ke1xXY2FsfXtcbWF0aGNhbHtXfX0KLS0tCgpgYGB7ciB4YXJpbmdhbkV4dHJhLWNsaXBib2FyZCwgZWNobyA9IEZBTFNFfQpodG1sdG9vbHM6OnRhZ0xpc3QoCiAgeGFyaW5nYW5FeHRyYTo6dXNlX2NsaXBib2FyZCgKICAgIGJ1dHRvbl90ZXh0ID0gIjxpIGNsYXNzPVwiZmEtc29saWQgZmEtY2xpcGJvYXJkXCIgc3R5bGU9XCJjb2xvcjogIzAwMDA4QlwiPjwvaT4iLAogICAgc3VjY2Vzc190ZXh0ID0gIjxpIGNsYXNzPVwiZmEgZmEtY2hlY2tcIiBzdHlsZT1cImNvbG9yOiAjOTBCRTZEXCI+PC9pPiIsCiAgICBlcnJvcl90ZXh0ID0gIjxpIGNsYXNzPVwiZmEgZmEtdGltZXMtY2lyY2xlXCIgc3R5bGU9XCJjb2xvcjogI0Y5NDE0NFwiPjwvaT4iCiAgKSwKICBybWFya2Rvd246Omh0bWxfZGVwZW5kZW5jeV9mb250X2F3ZXNvbWUoKQopCmBgYAoKCmBgYHtjc3MsIGVjaG8gPSBGQUxTRX0KYm9keSAubWFpbi1jb250YWluZXIgewogIG1heC13aWR0aDogMTAwJSAhaW1wb3J0YW50OwogIHdpZHRoOiAxMDAlICFpbXBvcnRhbnQ7Cn0KYm9keSB7CiAgbWF4LXdpZHRoOiAxMDAlICFpbXBvcnRhbnQ7Cn0KCmJvZHksIHRkIHsKICAgZm9udC1zaXplOiAxNnB4Owp9CmNvZGUucnsKICBmb250LXNpemU6IDE0cHg7Cn0KcHJlIHsKICBmb250LXNpemU6IDE0cHgKfQouY3VzdG9tLWJveCB7CiAgYmFja2dyb3VuZC1jb2xvcjogI2Y1ZjdmYTsgLyogTGlnaHQgZ3JleS1ibHVlIGJhY2tncm91bmQgKi8KICBib3JkZXItY29sb3I6ICNlMWU4ZWQ7IC8qIExpZ2h0IGJvcmRlciBjb2xvciAqLwogIGNvbG9yOiAjMmMzZTUwOyAvKiBEYXJrIHRleHQgY29sb3IgKi8KICBwYWRkaW5nOiAxNXB4OyAvKiBQYWRkaW5nIGluc2lkZSB0aGUgYm94ICovCiAgYm9yZGVyLXJhZGl1czogNXB4OyAvKiBSb3VuZGVkIGNvcm5lcnMgKi8KICBtYXJnaW4tYm90dG9tOiAyMHB4OyAvKiBTcGFjaW5nIGJlbG93IHRoZSBib3ggKi8KfQouY2FwdGlvbiB7CiAgbWFyZ2luOiBhdXRvOwogIHRleHQtYWxpZ246IGNlbnRlcjsKICBtYXJnaW4tYm90dG9tOiAyMHB4OyAvKiBTcGFjaW5nIGJlbG93IHRoZSBib3ggKi8KfQpgYGAKCgpMZXQgdXMgc2V0IHNvbWUgZ2xvYmFsIG9wdGlvbnMgZm9yIGFsbCBjb2RlIGNodW5rcyBpbiB0aGlzIGRvY3VtZW50LgoKCmBgYHtyfQojIFNldCBzZWVkIGZvciByZXByb2R1Y2liaWxpdHkKc2V0LnNlZWQoMTk4MikgCiMgU2V0IGdsb2JhbCBvcHRpb25zIGZvciBhbGwgY29kZSBjaHVua3MKa25pdHI6Om9wdHNfY2h1bmskc2V0KAogICMgRGlzYWJsZSBtZXNzYWdlcyBwcmludGVkIGJ5IFIgY29kZSBjaHVua3MKICBtZXNzYWdlID0gVFJVRSwgICAgCiAgIyBEaXNhYmxlIHdhcm5pbmdzIHByaW50ZWQgYnkgUiBjb2RlIGNodW5rcwogIHdhcm5pbmcgPSBUUlVFLCAgICAKICAjIFNob3cgUiBjb2RlIHdpdGhpbiBjb2RlIGNodW5rcyBpbiBvdXRwdXQKICBlY2hvID0gVFJVRSwgICAgICAgIAogICMgSW5jbHVkZSBib3RoIFIgY29kZSBhbmQgaXRzIHJlc3VsdHMgaW4gb3V0cHV0CiAgaW5jbHVkZSA9IFRSVUUsICAgICAKICAjIEV2YWx1YXRlIFIgY29kZSBjaHVua3MKICBldmFsID0gVFJVRSwgICAgICAgCiAgIyBFbmFibGUgY2FjaGluZyBvZiBSIGNvZGUgY2h1bmtzIGZvciBmYXN0ZXIgcmVuZGVyaW5nCiAgY2FjaGUgPSBGQUxTRSwgICAgICAKICAjIEFsaWduIGZpZ3VyZXMgaW4gdGhlIGNlbnRlciBvZiB0aGUgb3V0cHV0CiAgZmlnLmFsaWduID0gImNlbnRlciIsCiAgIyBFbmFibGUgcmV0aW5hIGRpc3BsYXkgZm9yIGhpZ2gtcmVzb2x1dGlvbiBmaWd1cmVzCiAgcmV0aW5hID0gMiwKICAjIFNob3cgZXJyb3JzIGluIHRoZSBvdXRwdXQgaW5zdGVhZCBvZiBzdG9wcGluZyByZW5kZXJpbmcKICBlcnJvciA9IFRSVUUsCiAgIyBEbyBub3QgY29sbGFwc2UgY29kZSBhbmQgb3V0cHV0IGludG8gYSBzaW5nbGUgYmxvY2sKICBjb2xsYXBzZSA9IEZBTFNFCikKIyBTdGFydCB0aGUgZmlndXJlIGNvdW50ZXIKZmlnX2NvdW50IDwtIDAKIyBEZWZpbmUgdGhlIGNhcHRpb25lciBmdW5jdGlvbgpjYXB0aW9uZXIgPC0gZnVuY3Rpb24oY2FwdGlvbikgewogIGZpZ19jb3VudCA8PC0gZmlnX2NvdW50ICsgMQogIHBhc3RlMCgiRmlndXJlICIsIGZpZ19jb3VudCwgIjogIiwgY2FwdGlvbikKfQojIERlZmluZSB0aGUgZnVuY3Rpb24gdG8gdHJ1bmNhdGUgYSBudW1iZXIgdG8gdHdvIGRlY2ltYWwgcGxhY2VzCnRydW5jYXRlX3RvX3R3byA8LSBmdW5jdGlvbih4KSB7CiAgZmxvb3IoeCAqIDEwMCkgLyAxMDAKfQoKbTF0YWJsZSA8LSByU1BERTo6Om0xdGFibGUKbTJ0YWJsZSA8LSByU1BERTo6Om0ydGFibGUKbTN0YWJsZSA8LSByU1BERTo6Om0zdGFibGUKbTR0YWJsZSA8LSByU1BERTo6Om00dGFibGUKYGBgCgoKYGBge3J9CiMgaW5zdGFsbC5wYWNrYWdlcygiSU5MQSIscmVwb3M9YyhnZXRPcHRpb24oInJlcG9zIiksSU5MQT0iaHR0cHM6Ly9pbmxhLnItaW5sYS1kb3dubG9hZC5vcmcvUi90ZXN0aW5nIiksIGRlcD1UUlVFKQojIGlubGEudXBncmFkZSh0ZXN0aW5nID0gVFJVRSkKIyByZW1vdGVzOjppbnN0YWxsX2dpdGh1YigiaW5sYWJydS1vcmcvaW5sYWJydSIsIHJlZiA9ICJkZXZlbCIpCiMgcmVtb3Rlczo6aW5zdGFsbF9naXRodWIoImRhdmlkYm9saW4vcnNwZGUiLCByZWYgPSAiZGV2ZWwiKQojIHJlbW90ZXM6Omluc3RhbGxfZ2l0aHViKCJkYXZpZGJvbGluL21ldHJpY2dyYXBoIiwgcmVmID0gImRldmVsIikKbGlicmFyeShJTkxBKQpsaWJyYXJ5KGlubGFicnUpCmxpYnJhcnkoclNQREUpCmxpYnJhcnkoTWV0cmljR3JhcGgpCmxpYnJhcnkoZ3JhdGVmdWwpCgpsaWJyYXJ5KHBsb3RseSkKYGBgCgoKIyBVdGlsaXRhcnkgZnVuY3Rpb25zCgpgYGB7cn0KIyBGb3IgZWFjaCBtIGFuZCBiZXRhLCB0aGlzIGZ1bmN0aW9uIHJldHVybnMgY19tL2Jfe20rMX0gYW5kIHRoZSByb290cyBvZiByYiBhbmQgcmMKbXkuZ2V0LnJvb3RzIDwtIGZ1bmN0aW9uKG9yZGVyLCBiZXRhKSB7CiAgbXQgPC0gZ2V0KHBhc3RlMCgibSIsIG9yZGVyLCAidGFibGUiKSkKICByYiA8LSByZXAoMCwgb3JkZXIgKyAxKQogIHJjIDwtIHJlcCgwLCBvcmRlcikKICBpZihvcmRlciA9PSAxKSB7CiAgICAgIHJjID0gYXBwcm94KG10JGJldGEsIG10W1twYXN0ZTAoInJjIildXSwgYmV0YSkkeQogICAgfSBlbHNlIHsKICAgICAgcmMgPSBzYXBwbHkoMTpvcmRlciwgZnVuY3Rpb24oaSkgewogICAgICAgIGFwcHJveChtdCRiZXRhLCBtdFtbcGFzdGUwKCJyYy4iLCBpKV1dLCBiZXRhKSR5CiAgICAgIH0pCiAgICB9CiAgICByYiA9IHNhcHBseSgxOihvcmRlcisxKSwgZnVuY3Rpb24oaSkgewogICAgICBhcHByb3gobXQkYmV0YSwgbXRbW3Bhc3RlMCgicmIuIiwgaSldXSwgeG91dCA9IGJldGEpJHkKICAgIH0pCiAgICBmYWN0b3IgPSBhcHByb3gobXQkYmV0YSwgbXQkZmFjdG9yLCB4b3V0ID0gYmV0YSkkeQogIHJldHVybihsaXN0KHJiID0gcmIsIHJjID0gcmMsIGZhY3RvciA9IGZhY3RvcikpCn0KCiMgR2l2ZW4gdGhlIHJvb3RzLCByZXR1cm4gdGhlIHBvbHlub21pYWwgY29lZmZpY2llbnRzIGluIGluY3JlYXNpbmcgb3JkZXIgbGlrZSBhK2J4K2N4XjIrLi4uCiMgU28gaWYgdGhlIG91dHB1dCBpcyBjKDYsIC01LCAxKSwgdGhhdCBtZWFucyA2IC0gNXggKyB4XjIKcG9seV9mcm9tX3Jvb3RzIDwtIGZ1bmN0aW9uKHJvb3RzKSB7CiAgY29lZiA8LSAxCiAgZm9yIChyIGluIHJvb3RzKSB7CiAgICBjb2VmIDwtIGNvbnZvbHZlKGNvZWYsIGMoMSwgLXIpLCB0eXBlID0gIm9wZW4iKQogIH0KICByZXR1cm4oY29lZikgIyByZXR1cm5zIGluIGluY3JlYXNpbmcgb3JkZXIgbGlrZSBhK2J4K2N4XjIrLi4uCn0KIAojIEdpdmVuIGEgc2V0IG9mIGNvbXBsZXggcm9vdHMsIHRoaXMgZnVuY3Rpb24gZ3JvdXBzIGNvbXBsZXggcm9vdHMgaW50byBjb25qdWdhdGUgcGFpcnMKZ2V0X2Nvbmp1Z2F0ZV9wYWlycyA8LSBmdW5jdGlvbihyb290cykgewogIHVzZWQgPC0gcmVwKEZBTFNFLCBsZW5ndGgocm9vdHMpKQogIHBhaXJzIDwtIGxpc3QoKQogIGZvciAoaSBpbiBzZXFfYWxvbmcocm9vdHMpKSB7CiAgICBpZiAoIXVzZWRbaV0pIHsKICAgICAgY29ual9yb290IDwtIENvbmoocm9vdHNbaV0pCiAgICAgICMgRmluZCBpdHMgY29uanVnYXRlICh3aXRoaW4gdG9sZXJhbmNlKQogICAgICBqIDwtIHdoaWNoKCF1c2VkICYgYWJzKFJlKHJvb3RzKSAtIFJlKGNvbmpfcm9vdCkpIDwgMWUtOCAmIGFicyhJbShyb290cykgLSBJbShjb25qX3Jvb3QpKSA8IDFlLTgpCiAgICAgIGogPC0galtqICE9IGldCiAgICAgIGlmIChsZW5ndGgoaikgPiAwKSB7CiAgICAgICAgdXNlZFtjKGksIGpbMV0pXSA8LSBUUlVFCiAgICAgICAgcGFpcnNbW2xlbmd0aChwYWlycykgKyAxXV0gPC0gYyhyb290c1tpXSwgcm9vdHNbalsxXV0pCiAgICAgIH0KICAgIH0KICB9CiAgcmV0dXJuKHBhaXJzKQp9CgojIFRoaXMgZnVuY3Rpb24gYnVpbGRzIHRoZSBwb2x5bm9taWFsIGNfbVxiX3ttKzF9IHByb2Rfe2k9MX1ebSAoMS1yX3sxaX14KSArIHRhdSBccHJvZF97aj0xfV57bSsxfSAoMS1yX3syan14KQojIEZvciBleGFtcGxlIGNvbXB1dGVfc3VtX3BvbHkoMCwgcmVwKDAsMiksIGMoMiwtMywxKSwxKSA9IGMoNiwtNywwLDEpIGJlY2F1c2UgKDEtMngpKDErM3gpKDEteCkgPSA2eF4zIC0gN3heMiArIDB4ICsgMQpjb21wdXRlX3N1bV9wb2x5IDwtIGZ1bmN0aW9uKGZhY3RvciwgcHJfcm9vdHMsIHBsX3Jvb3RzLCBjdGUpIHsKICBwcl9jb2VmIDwtIGMoMCwgcG9seV9mcm9tX3Jvb3RzKHByX3Jvb3RzKSkgIyBpbiBkZWNyZWFzaW5nIG9yZGVyIGxpa2UgeF5uK2J4XihuLTEpK2N4XihuLTIpKy4uLgogIHBsX2NvZWYgPC0gcG9seV9mcm9tX3Jvb3RzKHBsX3Jvb3RzKSAjIGluIGRlY3JlYXNpbmcgb3JkZXIgbGlrZSB4Xm4rYnheKG4tMSkrY3heKG4tMikrLi4uCiAgcHJfcGx1c19wbF9jb2VmIDwtIGZhY3RvciAqIHByX2NvZWYgKyBjdGUgKiBwbF9jb2VmCiAgcmV0dXJuKHByX3BsdXNfcGxfY29lZikgIyBpbiBkZWNyZWFzaW5nIG9yZGVyIGxpa2UgeF5uK2J4XihuLTEpK2N4XihuLTIpKy4uLgp9Cgpjb21wdXRlX3JlYWxfcm9vdHNfYW5kX2NvbXBsZXhfY29lZiA8LSBmdW5jdGlvbihwcl9wbHVzX3BsX2NvZWYpIHsKICBwcl9wbHVzX3BsX3Jvb3RzIDwtIHBvbHlyb290KHJldihwcl9wbHVzX3BsX2NvZWYpKSAjIHBvbHlyb290IGV4cGVjdHMgdGhlIGNvZWZmaWNpZW50cyBpbiBpbmNyZWFzaW5nIG9yZGVyCiAgcmVhbF9yb290cyA8LSBSZShwcl9wbHVzX3BsX3Jvb3RzW2FicyhJbShwcl9wbHVzX3BsX3Jvb3RzKSkgPCAxZS04XSkKICBjb21wbGV4X3Jvb3RzIDwtIHByX3BsdXNfcGxfcm9vdHNbYWJzKEltKHByX3BsdXNfcGxfcm9vdHMpKSA+PSAxZS04XQogIGNvbXBsZXhfcG9seV9jb2VmcyA8LSBsaXN0KCkKICBpZiAobGVuZ3RoKGNvbXBsZXhfcm9vdHMpID4gMCkgewogICAgcGFpcnMgPC0gZ2V0X2Nvbmp1Z2F0ZV9wYWlycyhjb21wbGV4X3Jvb3RzKQogICAgZm9yIChwYWlyIGluIHBhaXJzKSB7CiAgICAgIGNvZWYgPC0gcmV2KFJlKHBvbHlfZnJvbV9yb290cyhjKHBhaXJbMV0sIHBhaXJbMl0pKSkpICAjIHJldHVybnMgaW4geF4yICsgYnggKyBjIG9yZGVyCiAgICAgIGNvbXBsZXhfcG9seV9jb2Vmc1tbbGVuZ3RoKGNvbXBsZXhfcG9seV9jb2VmcykgKyAxXV0gPC0gcm91bmQoY29lZiwgMTIpCiAgICB9CiAgfSBlbHNlIHsKICAgIGNvbXBsZXhfcG9seV9jb2VmcyA8LSBsaXN0KCkgICMgTm8gY29tcGxleCByb290cwogIH0KICByZXR1cm4obGlzdChuZXdfZmFjdG9yID0gcHJfcGx1c19wbF9jb2VmWzFdLCAKICAgICAgICAgICAgICBwcl9wbHVzX3BsX3Jvb3RzID0gcHJfcGx1c19wbF9yb290cywgCiAgICAgICAgICAgICAgcmVhbF9yb290cyA9IHJlYWxfcm9vdHMsIAogICAgICAgICAgICAgIGNvbXBsZXhfcG9seV9jb2VmcyA9IGNvbXBsZXhfcG9seV9jb2VmcykpCn0KCm15LmZyYWN0aW9uYWwub3BlcmF0b3JzIDwtIGZ1bmN0aW9uKEwsICMga2FwcGFeMkMgKyBHCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGJldGEsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIEMsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHNjYWxlLmZhY3RvciwgIyBrYXBwYV4yCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIG0gPSAxLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aW1lX3N0ZXApIHsKCiAgQyA8LSBNYXRyaXg6OkRpYWdvbmFsKGRpbShDKVsxXSwgcm93U3VtcyhDKSkgIyBsdW1wZWQKICBDaSA8LSBNYXRyaXg6OkRpYWdvbmFsKGRpbShDKVsxXSwgMSAvIHJvd1N1bXMoQykpICMgbHVtcGVkIAogIEkgPC0gTWF0cml4OjpEaWFnb25hbChkaW0oQylbMV0pCiAgTCA8LSBMIC8gc2NhbGUuZmFjdG9yICMgQyArIEcva2FwcGFeMgogIExDaSA8LSBMICUqJSBDaQogIHJvb3RzIDwtIG15LmdldC5yb290cyhtLCBiZXRhKQogIFBsLnJvb3RzIDwtIHJvb3RzJHJiCiAgUHIucm9vdHMgPC0gcm9vdHMkcmMKICBmYWN0b3IgPC0gcm9vdHMkZmFjdG9yCiAgbmV3X2ZhY3Rvcl9hbmRfcm9vdHMgPC0gY29tcHV0ZV9yZWFsX3Jvb3RzX2FuZF9jb21wbGV4X2NvZWYoY29tcHV0ZV9zdW1fcG9seShmYWN0b3IsIFByLnJvb3RzLCBQbC5yb290cywgdGltZV9zdGVwKSkKICAKICBuZXdfcmVhbF9yb290cyA8LSBuZXdfZmFjdG9yX2FuZF9yb290cyRyZWFsX3Jvb3RzCiAgbmV3X2ZhY3RvciA8LSBuZXdfZmFjdG9yX2FuZF9yb290cyRuZXdfZmFjdG9yCiAgY29tcGxleF9wb2x5X2NvZWZzIDwtIG5ld19mYWN0b3JfYW5kX3Jvb3RzJGNvbXBsZXhfcG9seV9jb2VmcwogIAogICMgRmlsbCBQcl9wbHVzX1BsLmZhY3RvcnMKICBQcl9wbHVzX1BsLmZhY3RvcnMgPC0gbGlzdCgpCiAgaWYgKGxlbmd0aChuZXdfcmVhbF9yb290cykgPj0gMSkgewogICAgZm9yIChpIGluIDE6bGVuZ3RoKG5ld19yZWFsX3Jvb3RzKSkgewogICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbW2ldXSA8LSBMQ2kgLSBuZXdfcmVhbF9yb290c1tpXSAqIEkKICAgIH0KICB9CiAgCiAgcHJpbnQobGVuZ3RoKFByX3BsdXNfUGwuZmFjdG9ycykpCiAgCiAgbGVuZ3RoX25vdyA8LSBsZW5ndGgoUHJfcGx1c19QbC5mYWN0b3JzKQogIGlmKGxlbmd0aChjb21wbGV4X3BvbHlfY29lZnMpID49IDEpIHsKICAgIGZvciAoaSBpbiAxOmxlbmd0aChjb21wbGV4X3BvbHlfY29lZnMpKSB7CiAgICB9CiAgICBQcl9wbHVzX1BsLmZhY3RvcnNbW2xlbmd0aF9ub3cgKyBpXV0gPC0gTENpICUqJSBMQ2kgKyBjb21wbGV4X3BvbHlfY29lZnNbW2ldXVsyXSAqIExDaSArIGNvbXBsZXhfcG9seV9jb2Vmc1tbaV1dWzNdICogSQogIH0KICAjIEZpbGwgUHIuZmFjdG9ycwogIFByLmZhY3RvcnMgPC0gbGlzdCgpCiAgUHIuZmFjdG9yc1tbMV1dIDwtIEkgLSBMQ2kgKiByb290cyRyY1sxXQogIGlmIChsZW5ndGgocm9vdHMkcmMpID4gMSkgewogICAgZm9yIChpIGluIDI6bGVuZ3RoKHJvb3RzJHJjKSkgewogICAgICBQci5mYWN0b3JzW1tpXV0gPC0gSSAtIExDaSAqIHJvb3RzJHJjW2ldCiAgICB9CiAgfQogIAogIHJldHVybihsaXN0KAogICAgQ2kgPSBDaSwgQyA9IEMsIExDaSA9IExDaSwgTCA9IEwsIG0gPSBtLCBiZXRhID0gYmV0YSwKICAgIGZhY3RvciA9IGZhY3RvciwKICAgIG5ld19mYWN0b3IgPSBuZXdfZmFjdG9yLAogICAgUHIuZmFjdG9ycyA9IFByLmZhY3RvcnMsCiAgICBQcl9wbHVzX1BsLmZhY3RvcnMgPSBQcl9wbHVzX1BsLmZhY3RvcnMKICApKQp9CgpteS5zb2x2ZXIgPC0gZnVuY3Rpb24ob2JqLCB2KXsKICBQci5mYWN0b3JzIDwtIG9iaiRQci5mYWN0b3JzCiAgUHJfcGx1c19QbC5mYWN0b3JzIDwtIG9iaiRQcl9wbHVzX1BsLmZhY3RvcnMKICBtIDwtIG9iaiRtCiAgQyA8LSBvYmokQwogIENpIDwtIG9iaiRDaQogIGZhY3RvciA8LSBvYmokZmFjdG9yCiAgbmV3X2ZhY3RvciA8LSBvYmokbmV3X2ZhY3RvcgogIGlmIChtPT0xKXsKICAgIHRlbXAgPC0gc29sdmUoUHJfcGx1c19QbC5mYWN0b3JzW1syXV0sIAogICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1sxXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbMV1dLCAKICAgICAgICAgICAgICAgICAgICBDJSoldikpCiAgfSBlbHNlIGlmIChtPT0yKXsKICAgIHRlbXAgPC0gc29sdmUoUHJfcGx1c19QbC5mYWN0b3JzW1szXV0sIAogICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1syXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbMl1dLCAKICAgICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1sxXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgICAgUHJfcGx1c19QbC5mYWN0b3JzW1sxXV0sIAogICAgICAgICAgICAgICAgICAgICAgQyUqJXYpKSkKICB9IGVsc2UgaWYgKG09PTMpewogICAgdGVtcCA8LSBzb2x2ZShQcl9wbHVzX1BsLmZhY3RvcnNbWzRdXSwgCiAgICAgICAgICAgICAgICAgIFByLmZhY3RvcnNbWzNdXSAlKiUgc29sdmUoCiAgICAgICAgICAgICAgICAgICAgUHJfcGx1c19QbC5mYWN0b3JzW1szXV0sIAogICAgICAgICAgICAgICAgICAgIFByLmZhY3RvcnNbWzJdXSAlKiUgc29sdmUoCiAgICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzJdXSwgCiAgICAgICAgICAgICAgICAgICAgICBQci5mYWN0b3JzW1sxXV0gJSolIHNvbHZlKAogICAgICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzFdXSwgCiAgICAgICAgICAgICAgICAgICAgICAgIEMlKiV2KSkpKQogIH0gZWxzZSBpZiAobT09NCl7CiAgICB0ZW1wIDwtIHNvbHZlKFByX3BsdXNfUGwuZmFjdG9yc1tbNV1dLCAKICAgICAgICAgICAgICAgICAgUHIuZmFjdG9yc1tbNF1dICUqJSBzb2x2ZSgKICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzRdXSwgCiAgICAgICAgICAgICAgICAgICAgUHIuZmFjdG9yc1tbM11dICUqJSBzb2x2ZSgKICAgICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbM11dLCAKICAgICAgICAgICAgICAgICAgICAgIFByLmZhY3RvcnNbWzJdXSAlKiUgc29sdmUoCiAgICAgICAgICAgICAgICAgICAgICAgIFByX3BsdXNfUGwuZmFjdG9yc1tbMl1dLCAKICAgICAgICAgICAgICAgICAgICAgICAgUHIuZmFjdG9yc1tbMV1dICUqJSBzb2x2ZSgKICAgICAgICAgICAgICAgICAgICAgICAgICBQcl9wbHVzX1BsLmZhY3RvcnNbWzFdXSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgQyUqJXYpKSkpKQogIH0KICByZXR1cm4oKGZhY3Rvci9uZXdfZmFjdG9yKSAqIENpICUqJSB0ZW1wKQp9CmBgYAoKCldlIHdhbnQgdG8gc29sdmUgdGhlIGZyYWN0aW9uYWwgZGlmZnVzaW9uIGVxdWF0aW9uClxiZWdpbntlcXVhdGlvbn0KXGxhYmVse2VxOm1haW5lcX0KICAgIFxwYXJ0aWFsX3QgdSsoXGthcHBhXjItXERlbHRhX1xHYW1tYSlee1xhbHBoYS8yfSB1PWYgXHRleHQgeyBvbiB9IFxHYW1tYSBcdGltZXMoMCwgVCksIFxxdWFkIHUoMCk9dV8wIFx0ZXh0IHsgb24gfSBcR2FtbWEsClxlbmR7ZXF1YXRpb259CndoZXJlICR1JCBzYXRpc2ZpZXMgdGhlIEtpcmNoaG9mZiB2ZXJ0ZXggY29uZGl0aW9ucwpcYmVnaW57ZXF1YXRpb259ClxsYWJlbHtlcTpLY29uZH0KICAgIFxsZWZ0XHtccGhpXGluIEMoXEdhbW1hKVw7XEJpZ3xcOyBcZm9yYWxsIHZcaW4gVjogXHN1bV97ZVxpblxtYXRoY2Fse0V9X3Z9XHBhcnRpYWxfZSBccGhpKHYpPTAgXHJpZ2h0XH0KXGVuZHtlcXVhdGlvbn0KCklmICRmPTAkLCB0aGVuIHRoZSBzb2x1dGlvbiBpcyBnaXZlbiBieQpcYmVnaW57ZXF1YXRpb259ClxsYWJlbHtlcTpzb2xfcmVwcmVudGF0aW9ufQogICAgICAgIHUocyx0KSA9IFxkaXNwbGF5c3R5bGVcc3VtX3tqXGluXG1hdGhiYntOfX1lXnstXGxhbWJkYV57XGFscGhhLzJ9X2p0fVxsZWZ0KHVfMCwgZV9qXHJpZ2h0KV97TF8yKFxHYW1tYSl9ZV9qKHMpLgpcZW5ke2VxdWF0aW9ufQoKYGBge3J9CiMgRnVuY3Rpb24gdG8gYnVpbGQgYSB0YWRwb2xlIGdyYXBoIGFuZCBjcmVhdGUgYSBtZXNoCmdldHNfZ3JhcGhfdGFkcG9sZSA8LSBmdW5jdGlvbihoKXsKICBlZGdlMSA8LSByYmluZChjKDAsMCksYygxLDApKQogIHRoZXRhIDwtIHNlcShmcm9tPS1waSx0bz1waSxsZW5ndGgub3V0ID0gMTAwKQogIGVkZ2UyIDwtIGNiaW5kKDErMS9waStjb3ModGhldGEpL3BpLHNpbih0aGV0YSkvcGkpCiAgZWRnZXMgPSBsaXN0KGVkZ2UxLCBlZGdlMikKICBncmFwaCA8LSBtZXRyaWNfZ3JhcGgkbmV3KGVkZ2VzID0gZWRnZXMpCiAgZ3JhcGgkYnVpbGRfbWVzaChoID0gaCkKICByZXR1cm4oZ3JhcGgpCn0KYGBgCgpMZXQgJFxHYW1tYV9UID0gKFxWY2FsLFxFY2FsKSQgY2hhcmFjdGVyaXplIHRoZSB0YWRwb2xlIGdyYXBoIHdpdGggJFxWY2FsID0gXHt2XzEsdl8yXH0kIGFuZCAkXEVjYWwgPSBce2VfMSxlXzJcfSQuIFRoZSBsZWZ0IGVkZ2UgJGVfMSQgaGFzIGxlbmd0aCAxIGFuZCB0aGUgY2lyY3VsYXIgZWRnZSAkZV8yJCBoYXMgbGVuZ3RoIDIuIEFzIGRpc2N1c3NlZCBiZWZvcmUsIGEgcG9pbnQgb24gJGVfMSQgaXMgcGFyYW1ldGVyaXplZCB2aWEgJHM9XGxlZnQoZV8xLCB0XHJpZ2h0KSQgZm9yICR0IFxpblswLDFdJCBhbmQgYSBwb2ludCBvbiAkZV8yJCB2aWEgJHM9XGxlZnQoZV8yLCB0XHJpZ2h0KSQgZm9yICR0XGluWzAsMl0kLiBPbmUgY2FuIHZlcmlmeSB0aGF0ICQtXERlbHRhX1xHYW1tYSQgaGFzIGVpZ2VudmFsdWVzICQwLFxsZWZ0XHsoaSBccGkgLyAyKV4yXHJpZ2h0XH1fe2kgXGluIFxtYXRoYmJ7Tn19JCBhbmQgJFxsZWZ0XHsoaSBccGkgLyAyKV4yXHJpZ2h0XH1fezIgaSBcaW4gXG1hdGhiYntOfX0kIHdpdGggY29ycmVzcG9uZGluZyBlaWdlbmZ1bmN0aW9ucyAkXHBoaV8wJCwgJFxsZWZ0XHtccGhpX2lccmlnaHRcfV97aSBcaW4gXG1hdGhiYntOfX0kLCBhbmQgJFxsZWZ0XHtccHNpX2lccmlnaHRcfV97MiBpIFxpbiBcbWF0aGJie059fSQgZ2l2ZW4gYnkgJFxwaGlfMChzKT0xIC8gXHNxcnR7M30kIGFuZCAKXGJlZ2lue2VxdWF0aW9uKn0KICAgIFxwaGlfaShzKT1DX3tccGhpLCBpfVxiZWdpbntjYXNlc30KICAgICAgICAtMiBcc2luIChcZnJhY3tpXHBpfXsyfSkgXGNvcyAoXGZyYWN7aSBccGkgdH17Mn0pLCAmIHMgXGluIGVfMSwgXFwKXHNpbiAoaSBccGkgdCAvIDIpLCAmIHMgXGluIGVfMiwKICAgIFxlbmR7Y2FzZXN9LApccXVhZCAKICAgIFxwc2lfaShzKT1cZnJhY3tcc3FydHszfX17XHNxcnR7Mn19IFxiZWdpbntjYXNlc30KICAgICgtMSlee2kgLyAyfSBcY29zIChcZnJhY3tpIFxwaSB0fXsyfSksICYgcyBcaW4gZV8xLCBcXApcY29zIChcZnJhY3tpIFxwaSB0fXsyfSksICYgcyBcaW4gZV8yLApcZW5ke2Nhc2VzfSwKXGVuZHtlcXVhdGlvbip9CndoZXJlICRDX3tccGhpLCBpfT0xJCBpZiAkaSQgaXMgZXZlbiBhbmQgJENfe1xwaGksIGl9PTEgLyBcc3FydHszfSQgb3RoZXJ3aXNlLiBNb3Jlb3ZlciwgdGhlc2UgZnVuY3Rpb25zIGZvcm0gYW4gb3J0aG9ub3JtYWwgYmFzaXMgZm9yICRMXzIoXEdhbW1hX1QpJC4KCmBgYHtyfQojIEZ1bmN0aW9uIHRvIGNvbXB1dGUgdGhlIGVpZ2VuZnVuY3Rpb25zIAp0YWRwb2xlLmVpZyA8LSBmdW5jdGlvbihrLGdyYXBoKXsKeDEgPC0gYygwLGdyYXBoJGdldF9lZGdlX2xlbmd0aHMoKVsxXSpncmFwaCRtZXNoJFB0RVtncmFwaCRtZXNoJFB0RVssMV09PTEsMl0pIAp4MiA8LSBjKDAsZ3JhcGgkZ2V0X2VkZ2VfbGVuZ3RocygpWzJdKmdyYXBoJG1lc2gkUHRFW2dyYXBoJG1lc2gkUHRFWywxXT09MiwyXSkgCgppZihrPT0wKXsgCiAgZi5lMSA8LSByZXAoMSxsZW5ndGgoeDEpKSAKICBmLmUyIDwtIHJlcCgxLGxlbmd0aCh4MikpIAogIGYxID0gYyhmLmUxWzFdLGYuZTJbMV0sZi5lMVstMV0sIGYuZTJbLTFdKSAKICBmID0gbGlzdChwaGk9ZjEvc3FydCgzKSkgCiAgCn0gZWxzZSB7CiAgZi5lMSA8LSAtMipzaW4ocGkqayoxLzIpKmNvcyhwaSprKngxLzIpIAogIGYuZTIgPC0gc2luKHBpKmsqeDIvMikgICAgICAgICAgICAgICAgICAKICAKICBmMSA9IGMoZi5lMVsxXSxmLmUyWzFdLGYuZTFbLTFdLCBmLmUyWy0xXSkgCiAgCiAgaWYoKGsgJSUgMik9PTEpeyAKICAgIGYgPSBsaXN0KHBoaT1mMS9zcXJ0KDMpKSAKICB9IGVsc2UgeyAKICAgIGYuZTEgPC0gKC0xKV57ay8yfSpjb3MocGkqayp4MS8yKQogICAgZi5lMiA8LSBjb3MocGkqayp4Mi8yKQogICAgZjIgPSBjKGYuZTFbMV0sZi5lMlsxXSxmLmUxWy0xXSxmLmUyWy0xXSkgCiAgICBmIDwtIGxpc3QocGhpPWYxLHBzaT1mMi9zcXJ0KDMvMikpCiAgfQp9CgpyZXR1cm4oZikKfQpgYGAKCkxldCAkXGFscGhhXGluKDAsMl0kIGFuZCAkVV9oXlx0YXUkIGRlbm90ZSB0aGUgc2VxdWVuY2Ugb2YgYXBwcm94aW1hdGlvbnMgb2YgdGhlIHNvbHV0aW9uIHRvIHRoZSB3ZWFrIGZvcm0gb2YgcHJvYmxlbSBcZXFyZWZ7ZXE6bWFpbmVxfSBhdCBlYWNoIHRpbWUgc3RlcCBvbiBhIG1lc2ggaW5kZXhlZCBieSAkaCQuIExldCAkej0wJCBhbmQgJFVeMF9oID0gUF9odV8wJC4gRm9yICRrPTAsXGRvdHMsIE4tMSQsICRVX2hee2srMX1caW4gVl9oJCBzb2x2ZXMgdGhlIGZvbGxvd2luZyBzY2hlbWUKXGJlZ2lue2FsaWdufQpcbGFiZWx7c3lzdGVtOmZ1bGx5X2Rpc2NyZXRlX3NjaGVtZX0KICAgICAgICBcbGFuZ2xlXGRlbHRhIFVfaF57aysxfSxccGhpXHJhbmdsZSArIFxtYXRoZnJha3thfShVX2hee2srMX0sXHBoaSkgPSBcbGFuZ2xlIGZee2srMX0sXHBoaVxyYW5nbGUgLFxxdWFkXGZvcmFsbFxwaGlcaW4gVl9oLApcZW5ke2FsaWdufQp3aGVyZSAkZl57aysxfSA9IFxkaXNwbGF5c3R5bGVcZGZyYWN7MX17XHRhdX1caW50X3t0X2t9Xnt0XntrKzF9fWYodClkdCQuClRoZSBzb2x1dGlvbiBjYW4gYmUgcmVwcmVzZW50ZWQgYXMgClxiZWdpbnthbGlnbip9CiAgICBVX2heayhzKSA9ICBcc3VtX3tqPTF9XntOX2h9dV9qXmtccHNpX2oocykKXGVuZHthbGlnbip9ClJlcGxhY2luZyB0aGlzIGludG8gXGVxcmVme3N5c3RlbTpmdWxseV9kaXNjcmV0ZV9zY2hlbWV9IHlpZWxkcyB0aGUgZm9sbG93aW5nIGxpbmVhciBzeXN0ZW0KXGJlZ2lue2FsaWduKn0KICAgIFxzdW1fe2o9MX1ee05faH11X2pee2srMX1bKFxwc2lfaixccHNpX2kpX3tMXzIoXEdhbW1hKX0rIFx0YXVcbWF0aGZyYWt7YX0oXHBzaV9qLFxwc2lfaSldID0gXHN1bV97aj0xfV57Tl9ofXVfal57a30oXHBzaV9qLFxwc2lfaSlfe0xfMihcR2FtbWEpfStcdGF1KCBmXntrKzF9LFxwc2lfaSlfe0xfMihcR2FtbWEpfQpcZW5ke2FsaWduKn0KZm9yICRpID0gMSxcZG90cywgTl9oJC4gSW4gbWF0cml4IG5vdGF0aW9uLApcYmVnaW57YWxpZ259ClxsYWJlbHtkaWZmX2VxX2Rpc2NyZXRlfQogICAgKEMrXHRhdSBMXntcYWxwaGEvMn0pVV57aysxfSA9IENVXmsrXHRhdSBGXntrKzF9LApcZW5ke2FsaWdufQp3aGVyZSAkQyQgaGFzIGVudHJpZXMgJENfe2lqfSA9IChccHNpX2osXHBzaV9pKV97TF8yKFxHYW1tYSl9JCwgJExee1xhbHBoYS8yfSQgaGFzIGVudHJpZXMgJFxtYXRoZnJha3thfShccHNpX2osXHBzaV9pKSQsICRVXmskIGhhcyBlbnRyaWVzICR1X2peayQsIGFuZCAkRl5rJCBoYXMgZW50cmllcyAkKCBmXntrfSxccHNpX2kpX3tMXzIoXEdhbW1hKX0kLiBGb3Igc2ltcGxpY2l0eSwgbGV0ICRmPTAkIGFuZCBzbyB3ZSBhcnJpdmUgYXQKClxiZWdpbntlcXVhdGlvbn0KKEMrXHRhdSBMXntcYWxwaGEvMn0pVV57aysxfSA9IENVXntrfQpcZW5ke2VxdWF0aW9ufQoKQXBwbHkgJExeey1cYWxwaGEvMn0kIHRvIGJvdGggc2lkZXMgdG8gZ2V0CgpcYmVnaW57ZXF1YXRpb259CihMXnstXGFscGhhLzJ9QytcdGF1IEkpVV57aysxfSA9IExeey1cYWxwaGEvMn1DVV57a30KXGVuZHtlcXVhdGlvbn0KCkFwcHJveGltYXRlICRMXnstXGFscGhhLzJ9JCBieSAkUF9yUF9sXnstMX0kIHRvIGFycml2ZSBhdAoKXGJlZ2lue2VxdWF0aW9ufQooUF9yUF9sXnstMX1DK1x0YXUgSSlVXntrKzF9ID0gUF9yUF9sXnstMX1DVV57a30KXGVuZHtlcXVhdGlvbn0KCkFjdHVhbGx5LCAkTF57LVxhbHBoYS8yfSQgc2hvdWxkIGJlIGFwcHJveGltYXRlZCBieSAkUF9sXnstVH1QX3JeVCQgKHRoaXMgaXMgbGlrZSB0aGlzIGJlY2F1c2Ugd2Ugd2FudCB0byB1c2UgJFBfciQgYW5kICRQX2wkIHRoZSB3YXkgdGhleSB3aGVyZSBjb25zdHJ1Y3RlZCkgYW5kIHNvIHdlIGNhbiB3cml0ZQoKXGJlZ2lue2VxdWF0aW9ufQooUF9sXnstVH1QX3JeVEMrXHRhdSBJKVVee2srMX0gPSBQX2xeey1UfVBfcl5UQ1Vee2t9ClxlbmR7ZXF1YXRpb259CgpGcm9tIHRoaXMsIHdlIGFycml2ZSBhdCB0aGUgc2NoZW1lCgpcYmVnaW57ZXF1YXRpb259CihQX3JeVEMrXHRhdSBQX1xlbGxeVClVXntrKzF9ID0gUF9yXlRDVV5rClxsYWJlbHtlcTpzY2hlbWV9ClxlbmR7ZXF1YXRpb259CndoZXJlIChoZXJlIHdlIGFyZSB1c2luZyAkUF9yJCBhbmQgJFBfbCQgdGhlIHdheSB0aGV5IHdlcmUgY29uc3RydWN0ZWQpClxiZWdpbntlcXVhdGlvbn0KUF9yID0gY19tXHByb2Rfe2k9MX1ebSAoSS1yX3sxaX1DXnstMX1MKVxxdWFkXHRleHR7YW5kfVxxdWFkIFBfXGVsbCA9IGJfe20rMX1DXHByb2Rfe2o9MX1ee20rMX0gKEktcl97Mmp9Q157LTF9TCkKXGVuZHtlcXVhdGlvbn0KT2JzZXJ2ZSB0aGF0IApcYmVnaW57ZXF1YXRpb259ClBfcl5UID0gY19tXHByb2Rfe2k9MX1ebSAoSS1yX3sxaX1MQ157LTF9KVxxdWFkXHRleHR7YW5kfVxxdWFkIFBfXGVsbF5UID0gYl97bSsxfVxwcm9kX3tqPTF9XnttKzF9IChJLXJfezJqfUxDXnstMX0pXGNkb3QgQwpcZW5ke2VxdWF0aW9ufQpzaW5jZSAkTCQgYW5kICRDXnstMX0kIGFyZSBzeW1tZXRyaWMgYW5kIHRoZSBmYWN0b3JzIGluIHRoZSBwcm9kdWN0IGNvbW11dGUgKGJlY2F1c2UgdGhlIGZhY3RvcnMgaW4gdGhlIHByb2R1Y3QgdGhhdCBtYWtlIHVwIGEgcG9seW5vbWlhbCBjb21tdXRlLCBmb3IgZXhhbXBsZSwgJCgxLWF4KSgxLWJ4KSA9ICgxLWJ4KSgxLWF4KSQpLiBSZXBsYWNpbmcgdGhlc2UgdHdvIGluIG91ciBzY2hlbWUgd2UgZ2V0ClxiZWdpbntlcXVhdGlvbn0KXGxlZnQoXGRmcmFje2NfbX17Yl97bSsxfX1ccHJvZF97aT0xfV5tIChJLXJfezFpfUxDXnstMX0pK1x0YXUgXHByb2Rfe2o9MX1ee20rMX0gKEktcl97Mmp9TENeey0xfSlccmlnaHQpQ1Vee2srMX0gPSBcZGZyYWN7Y19tfXtiX3ttKzF9fVxwcm9kX3tpPTF9Xm0gKEktcl97MWl9TENeey0xfSlcY2RvdCBDVV5rClxlbmR7ZXF1YXRpb259CgpXZSBub3cgbmVlZCB0byBjb21wdXRlIHRoZSByb290cyBvZiB0aGUgcG9seW5vbWlhbApcYmVnaW57ZXF1YXRpb259ClIoeCkgPSBcZGZyYWN7Y19tfXtiX3ttKzF9fSBccHJvZF97aT0xfV5tICgxLXJfezFpfXgpK1x0YXUgXHByb2Rfe2o9MX1ee20rMX0gKDEtcl97Mmp9eCkgClxlbmR7ZXF1YXRpb259CiRSJCBoYXMgbGVhZGluZyBjb2VmZmljaWVudCAkXGdhbW1hPVx0YXUgKC0xKV57bSsxfSBccHJvZF97aj0xfV57bSsxfXJfezJqfSQgYW5kICRtKzEkIHJvb3RzICRSX2skIGZvciAkaz0xLFxsZG90cyxtKzEkLiBUaGF0IGlzLApcYmVnaW57ZXF1YXRpb259ClIoeCkgPSBcdGF1KC0xKV57bSsxfSBccHJvZF97aj0xfV57bSsxfXJfezJqfVxwcm9kX3trPTF9XnttKzF9ICh4LVJfayk9IFxnYW1tYVxwcm9kX3trPTF9XnttKzF9ICh4LVJfaykKXGVuZHtlcXVhdGlvbn0KV2UgY2FuIHRoZW4gd3JpdGUgb3VyIHNjaGVtZSBhcyAKXGJlZ2lue2VxdWF0aW9ufQpcbGVmdChcZ2FtbWFccHJvZF97az0xfV57bSsxfSAoTENeey0xfS1SX2tJKVxyaWdodClDVV57aysxfSA9IFxkZnJhY3tjX219e2Jfe20rMX19XHByb2Rfe2k9MX1ebSAoSS1yX3sxaX1MQ157LTF9KVxjZG90IENVXmsKXGVuZHtlcXVhdGlvbn0KVGhhdCBpcywKXGJlZ2lue2VxdWF0aW9ufQpVXntrKzF9ID0gXGRmcmFje2NfbX17Yl97bSsxfX1cZGZyYWN7MX17XGdhbW1hfUNeey0xfVxwcm9kX3trPTF9XnttKzF9IChMQ157LTF9LVJfa0kpXnstMX1ccHJvZF97aT0xfV5tIChJLXJfezFpfUxDXnstMX0pXGNkb3QgQ1VeawpcZW5ke2VxdWF0aW9ufQoKClJvb3RzICRyX3sxaX0kIGFuZCAkcl97Mmp9JCBhcmUgcmVhbCBidXQgJFJfayQgY2FuIGJlIGNvbXBsZXguIEluIHRoYXQgY2FzZSwgd2UgY29uc2lkZXIgdGhlIHF1YWRyYXRpYyB0ZXJtcy4gRm9yIGV4YW1wbGUsIGlmICRtPTMkLCBhbmQgJFIkIGhhcyBvbmUgcmVhbCByb290ICRyJCBhbmQgdHdvIGNvbXBsZXggY29uanVnYXRlIHJvb3RzICR6JCBhbmQgJFxiYXJ7en0kIHN1Y2ggdGhhdCAkKHgteikoeC1cYmFye3p9KSA9IHheMitheCtiJCwgd2Ugd2lsbCB3b3JrIHdpdGggClxiZWdpbntlcXVhdGlvbn0KUih4KSA9IFxnYW1tYSh4XjIrYXgrYikoeC1yKQpcZW5ke2VxdWF0aW9ufQppbnN0ZWFkIG9mIApcYmVnaW57ZXF1YXRpb259ClIoeCkgPSBcZ2FtbWEoeC16KSh4LVxiYXJ7en0pKHgtcikKXGVuZHtlcXVhdGlvbn0KCmBgYHtyfQpoIDwtIDAuMDEKZ3JhcGggPC0gZ2V0c19ncmFwaF90YWRwb2xlKGggPSBoKQpUX2ZpbmFsIDwtIDIKdGltZV9zdGVwIDwtIDAuMDEKdGltZV9zZXEgPC0gc2VxKDAsIFRfZmluYWwsIGJ5ID0gdGltZV9zdGVwKQojIENvbXB1dGUgdGhlIEZFTSBtYXRyaWNlcwpncmFwaCRjb21wdXRlX2ZlbSgpCkcgPC0gZ3JhcGgkbWVzaCRHCkMgPC0gZ3JhcGgkbWVzaCRDCkkgPC0gTWF0cml4OjpEaWFnb25hbChucm93KEMpKQp4IDwtIGdyYXBoJG1lc2gkVlssIDFdCnkgPC0gZ3JhcGgkbWVzaCRWWywgMl0KZWRnZV9udW1iZXIgPC0gZ3JhcGgkbWVzaCRWdEVbLCAxXQpwb3MgPC0gc3VtKGVkZ2VfbnVtYmVyID09IDEpKzEKb3JkZXJfdG9fcGxvdCA8LSBmdW5jdGlvbih2KXJldHVybihjKHZbMV0sIHZbMzpwb3NdLCB2WzJdLCB2Wyhwb3MrMSk6bGVuZ3RoKHYpXSwgdlsyXSkpCndlaWdodHMgPC0gZ3JhcGgkbWVzaCR3ZWlnaHRzCmBgYAoKCmBgYHtyfQprYXBwYSA8LSAxCmFscGhhIDwtIDEuOTkgIyBmcm9tIDAuNSB0byAyCm0gPSA0CmJldGEgPC0gYWxwaGEvMgpMIDwtIGthcHBhXjIqQyArIEcKYGBgCgojIFRydWUgc29sdXRpb24KCmBgYHtyfQojIEluaXRpYWwgY29uZGl0aW9uClVfMCA8LSAxMCpleHAoLSgoeC0xKV4yICsgKHkpXjIpKQoKVV90cnVlIDwtIG1hdHJpeChOQSwgbnJvdyA9IG5yb3coQyksIG5jb2wgPSBsZW5ndGgodGltZV9zZXEpKQpVX3RydWVbLCAxXSA8LSBVXzAKbl9maW5pdGUgPC0gMTAwCgoKZm9yIChrIGluIDE6KGxlbmd0aCh0aW1lX3NlcSkgLSAxKSkgewogIGF1eF9rIDwtIHJlcCgwLCBucm93KEMpKQogIGZvciAoaiBpbiAwOm5fZmluaXRlKSB7CiAgICBkZWNheV9qIDwtIGV4cCgtdGltZV9zZXFbaysxXSooa2FwcGFeMiArIChqKnBpLzIpXjIpXihhbHBoYS8yKSkKICAgIGVfaiA8LSB0YWRwb2xlLmVpZyhqLGdyYXBoKSRwaGkKICAgIGF1eF9rIDwtIGF1eF9rICsgZGVjYXlfaipzdW0oVV8wKmVfaip3ZWlnaHRzKSplX2oKICAgIGlmIChqPjAgJiYgKGogJSUgMiA9PSAwKSkgewogICAgICBlX2ogPC0gdGFkcG9sZS5laWcoaixncmFwaCkkcHNpCiAgICAgIGF1eF9rIDwtIGF1eF9rICsgZGVjYXlfaipzdW0oVV8wKmVfaip3ZWlnaHRzKSplX2oKICAgICAgfQogICAgfQogIFVfdHJ1ZVssIGsgKyAxXSA8LSBhdXhfawp9CmBgYAoKYGBge3J9Cm15X29wIDwtIG15LmZyYWN0aW9uYWwub3BlcmF0b3JzKEwsIGJldGEsIEMsIHNjYWxlLmZhY3RvciA9IGthcHBhXjIsIG0gPSBtLCB0aW1lX3N0ZXApCgpVX2FwcHJveDEgPC0gbWF0cml4KE5BLCBucm93ID0gbnJvdyhDKSwgbmNvbCA9IGxlbmd0aCh0aW1lX3NlcSkpClVfYXBwcm94MVssIDFdIDwtIFVfMAoKIyBUaW1lLXN0ZXBwaW5nIGxvb3AKZm9yIChrIGluIDE6KGxlbmd0aCh0aW1lX3NlcSkgLSAxKSkgewogIFVfYXBwcm94MVssIGsgKyAxXSA8LSBhcy5tYXRyaXgobXkuc29sdmVyKG15X29wLCBVX2FwcHJveDFbLCBrXSkpCn0KYGBgCgpgYGAKe3J9Cm9wIDwtIGZyYWN0aW9uYWwub3BlcmF0b3JzKEwsIGJldGEsIEMsIHNjYWxlLmZhY3RvciA9IGthcHBhXjIsIG0gPSBtKQpQbCA8LSBvcCRQbApQciA8LSBvcCRQcgpDaSA8LSBvcCRDaQpDIDwtIG9wJEMKYGBgCgoKYGBgCntyfQpMSFMgPC0gdChQciklKiUgQyArIHRpbWVfc3RlcCAqIHQoUGwpCiMgSW5pdGlhbGl6ZSBVIG1hdHJpeCB0byBzdG9yZSBzb2x1dGlvbiBhdCBlYWNoIHRpbWUgc3RlcApVX2FwcHJveDIgPC0gbWF0cml4KE5BLCBucm93ID0gbnJvdyhDKSwgbmNvbCA9IGxlbmd0aCh0aW1lX3NlcSkpClVfYXBwcm94MlssIDFdIDwtIFVfMAoKIyBUaW1lLXN0ZXBwaW5nIGxvb3AKZm9yIChrIGluIDE6KGxlbmd0aCh0aW1lX3NlcSkgLSAxKSkgewogICMgQ29tcHV0ZSB0aGUgcmlnaHQtaGFuZCBzaWRlIGZvciB0aGUgc2Vjb25kIGVxdWF0aW9uCiAgUkhTIDwtIHQoUHIpJSolIEMgJSolIFVfYXBwcm94MlssIGtdCiAgVV9hcHByb3gyWywgayArIDFdIDwtIGFzLm1hdHJpeChzb2x2ZShMSFMsIFJIUykpCn0KYGBgCgpgYGAKe3J9ClByLmFwcGx5Lm11bHQgPC0gZnVuY3Rpb24odil7cmV0dXJuKFByLm11bHQob3AsIHYpKX0KUGwuYXBwbHkuc29sdmUgPC0gZnVuY3Rpb24odil7cmV0dXJuKFBsLnNvbHZlKG9wLCB2KSl9ClBsaUMgPC0gYXBwbHkoQywgMiwgUGwuYXBwbHkuc29sdmUpICMgUGxDXi0xCiMgUHJlY29tcHV0ZSB0aGUgTEhTMSBtYXRyaXgKYXV4IDwtIGFwcGx5KFBsaUMsIDIsIFByLmFwcGx5Lm11bHQpICNQclBsXi0xQwpMSFMgPC0gYXV4ICsgdGltZV9zdGVwICogTWF0cml4OjpEaWFnb25hbChucm93KEMpKSAKCiMgSW5pdGlhbGl6ZSBVIG1hdHJpeCB0byBzdG9yZSBzb2x1dGlvbiBhdCBlYWNoIHRpbWUgc3RlcApVX2FwcHJveDIgPC0gbWF0cml4KE5BLCBucm93ID0gbnJvdyhDKSwgbmNvbCA9IGxlbmd0aCh0aW1lX3NlcSkpClVfYXBwcm94MlssIDFdIDwtIFVfMAoKIyBUaW1lLXN0ZXBwaW5nIGxvb3AKZm9yIChrIGluIDE6KGxlbmd0aCh0aW1lX3NlcSkgLSAxKSkgewogICMgQ29tcHV0ZSB0aGUgcmlnaHQtaGFuZCBzaWRlIGZvciB0aGUgc2Vjb25kIGVxdWF0aW9uCiAgUkhTIDwtIGF1eCAlKiUgVV9hcHByb3gyWywga10KICBVX2FwcHJveDJbLCBrICsgMV0gPC0gYXMubWF0cml4KHNvbHZlKExIUywgUkhTKSkKfQpgYGAKCiMgUGxvdAoKCmBgYHtyfQp4IDwtIG9yZGVyX3RvX3Bsb3QoeCkKeSA8LSBvcmRlcl90b19wbG90KHkpCgpVX2FwcHJveDIgPC0gVV9hcHByb3gxCgptYXhfZXJyb3JfYXRfZWFjaF90aW1lMSA8LSBhcHBseSgoVV90cnVlIC0gVV9hcHByb3gxKV4yLCAyLCBtZWFuKQptYXhfZXJyb3JfYXRfZWFjaF90aW1lMiA8LSBhcHBseSgoVV90cnVlIC0gVV9hcHByb3gyKV4yLCAyLCBtZWFuKQptYXhfZXJyb3JfYmV0d2Vlbl9ib3RoX2FwcHJveCA8LSBhcHBseSgoVV9hcHByb3gxIC0gVV9hcHByb3gyKV4yLCAyLCBtZWFuKQoKVV90cnVlIDwtIGFwcGx5KFVfdHJ1ZSwgMiwgb3JkZXJfdG9fcGxvdCkKVV9hcHByb3gxIDwtIGFwcGx5KFVfYXBwcm94MSwgMiwgb3JkZXJfdG9fcGxvdCkKVV9hcHByb3gyIDwtIGFwcGx5KFVfYXBwcm94MiwgMiwgb3JkZXJfdG9fcGxvdCkKCiMgQ3JlYXRlIGludGVyYWN0aXZlIHBsb3QKZmlnIDwtIHBsb3RfbHkoKQoKIyBBZGQgZmlyc3QgbGluZSAobWF4X2Vycm9yX2F0X2VhY2hfdGltZTEpCmZpZyA8LSBmaWcgJT4lIGFkZF90cmFjZSgKICB4ID0gfnRpbWVfc2VxLCB5ID0gfm1heF9lcnJvcl9hdF9lYWNoX3RpbWUxLCB0eXBlID0gJ3NjYXR0ZXInLCBtb2RlID0gJ2xpbmVzK21hcmtlcnMnLAogIGxpbmUgPSBsaXN0KGNvbG9yID0gJ3JlZCcsIHdpZHRoID0gMiksCiAgbWFya2VyID0gbGlzdChjb2xvciA9ICdyZWQnLCBzaXplID0gNCksCiAgbmFtZSA9ICJFcnJvciAgVHJ1ZSBhbmQgQXBwcm94IDEiCikKCiMgQWRkIHNlY29uZCBsaW5lIChtYXhfZXJyb3JfYXRfZWFjaF90aW1lMikKZmlnIDwtIGZpZyAlPiUgYWRkX3RyYWNlKAogIHggPSB+dGltZV9zZXEsIHkgPSB+bWF4X2Vycm9yX2F0X2VhY2hfdGltZTIsIHR5cGUgPSAnc2NhdHRlcicsIG1vZGUgPSAnbGluZXMrbWFya2VycycsCiAgbGluZSA9IGxpc3QoY29sb3IgPSAnYmx1ZScsIHdpZHRoID0gMiksCiAgbWFya2VyID0gbGlzdChjb2xvciA9ICdibHVlJywgc2l6ZSA9IDQpLAogIG5hbWUgPSAiRXJyb3IgVHJ1ZSBhbmQgQXBwcm94IDIiCikKCiMgQWRkIHRoaXJkIGxpbmUgKG1heF9lcnJvcl9iZXR3ZWVuX2JvdGhfYXBwcm94KQoKZmlnIDwtIGZpZyAlPiUgYWRkX3RyYWNlKAogIHggPSB+dGltZV9zZXEsIHkgPSB+bWF4X2Vycm9yX2JldHdlZW5fYm90aF9hcHByb3gsIHR5cGUgPSAnc2NhdHRlcicsIG1vZGUgPSAnbGluZXMrbWFya2VycycsCiAgbGluZSA9IGxpc3QoY29sb3IgPSAnb3JhbmdlJywgd2lkdGggPSAyKSwKICBtYXJrZXIgPSBsaXN0KGNvbG9yID0gJ29yYW5nZScsIHNpemUgPSA0KSwKICBuYW1lID0gIkVycm9yIEJldHdlZW4gQXBwcm94aW1hdGlvbnMiCikKCiMgTGF5b3V0CmZpZyA8LSBmaWcgJT4lIGxheW91dCgKICB0aXRsZSA9ICJFcnJvciBhdCBFYWNoIFRpbWUgU3RlcCIsCiAgeGF4aXMgPSBsaXN0KHRpdGxlID0gIlRpbWUiKSwKICB5YXhpcyA9IGxpc3QodGl0bGUgPSAiRXJyb3IiKSwKICBsZWdlbmQgPSBsaXN0KHggPSAwLjEsIHkgPSAwLjkpCikKCgpwbG90X2RhdGEgPC0gZGF0YS5mcmFtZSgKICB4ID0gcmVwKHgsIHRpbWVzID0gbmNvbChVX3RydWUpKSwKICB5ID0gcmVwKHksIHRpbWVzID0gbmNvbChVX3RydWUpKSwKICB6X3RydWUgPSBhcy52ZWN0b3IoVV90cnVlKSwKICB6X2FwcHJveDEgPSBhcy52ZWN0b3IoVV9hcHByb3gxKSwKICB6X2FwcHJveDIgPSBhcy52ZWN0b3IoVV9hcHByb3gyKSwKICBmcmFtZSA9IHJlcCh0aW1lX3NlcSwgZWFjaCA9IGxlbmd0aCh4KSkKKQoKIyBDb21wdXRlIGF4aXMgbGltaXRzCnhfcmFuZ2UgPC0gcmFuZ2UoeCkKeV9yYW5nZSA8LSByYW5nZSh5KQp6X3JhbmdlIDwtIHJhbmdlKGMoVV90cnVlLCBVX2FwcHJveDEsIFVfYXBwcm94MikpCgojIEluaXRpYWwgcGxvdCBzZXR1cCAoZmlyc3QgZnJhbWUgb25seSkKcCA8LSBwbG90X2x5KHBsb3RfZGF0YSwgZnJhbWUgPSB+ZnJhbWUpICU+JQogIGFkZF90cmFjZSgKICAgIHggPSB+eCwgeSA9IH55LCB6ID0gfnpfdHJ1ZSwKICAgIHR5cGUgPSAic2NhdHRlcjNkIiwgbW9kZSA9ICJsaW5lcyIsCiAgICBuYW1lID0gIlRydWUiLAogICAgbGluZSA9IGxpc3QoY29sb3IgPSAiZ3JlZW4iLCB3aWR0aCA9IDIpCiAgKSAlPiUKICBhZGRfdHJhY2UoCiAgICB4ID0gfngsIHkgPSB+eSwgeiA9IH56X2FwcHJveDEsCiAgICB0eXBlID0gInNjYXR0ZXIzZCIsIG1vZGUgPSAibGluZXMiLAogICAgbmFtZSA9ICJBcHByb3ggMSIsCiAgICBsaW5lID0gbGlzdChjb2xvciA9ICJyZWQiLCB3aWR0aCA9IDIpCiAgKSAlPiUKICBhZGRfdHJhY2UoCiAgICB4ID0gfngsIHkgPSB+eSwgeiA9IH56X2FwcHJveDIsCiAgICB0eXBlID0gInNjYXR0ZXIzZCIsIG1vZGUgPSAibGluZXMiLAogICAgbmFtZSA9ICJBcHByb3ggMiIsCiAgICBsaW5lID0gbGlzdChjb2xvciA9ICJibHVlIiwgd2lkdGggPSAyKQogICkgJT4lCiAgbGF5b3V0KAogICAgc2NlbmUgPSBsaXN0KAogICAgICB4YXhpcyA9IGxpc3QodGl0bGUgPSAieCIsIHJhbmdlID0geF9yYW5nZSksCiAgICAgIHlheGlzID0gbGlzdCh0aXRsZSA9ICJ5IiwgcmFuZ2UgPSB5X3JhbmdlKSwKICAgICAgemF4aXMgPSBsaXN0KHRpdGxlID0gIlZhbHVlIiwgcmFuZ2UgPSB6X3JhbmdlKSwKICAgICAgYXNwZWN0cmF0aW8gPSBsaXN0KHggPSAyLjQsIHkgPSAxLjIsIHogPSAxLjIpLAogICAgICAgICAgIGNhbWVyYSA9IGxpc3QoCiAgICAgIGV5ZSA9IGxpc3QoeCA9IDEuNSwgeSA9IDEuNSwgeiA9IDEpLCAgIyBBZGp1c3QgdGhlIHZpZXdwb2ludAogICAgICBjZW50ZXIgPSBsaXN0KHggPSAwLCB5ID0gMCwgeiA9IDApKQogICAgKSwKICAgIHVwZGF0ZW1lbnVzID0gbGlzdCgKICAgICAgbGlzdCgKICAgICAgICB0eXBlID0gImJ1dHRvbnMiLCBzaG93YWN0aXZlID0gRkFMU0UsCiAgICAgICAgYnV0dG9ucyA9IGxpc3QoCiAgICAgICAgICBsaXN0KGxhYmVsID0gIlBsYXkiLCBtZXRob2QgPSAiYW5pbWF0ZSIsCiAgICAgICAgICAgICAgIGFyZ3MgPSBsaXN0KE5VTEwsIGxpc3QoZnJhbWUgPSBsaXN0KGR1cmF0aW9uID0gMTAwLCByZWRyYXcgPSBUUlVFKSwgZnJvbWN1cnJlbnQgPSBUUlVFKSkpLAogICAgICAgICAgbGlzdChsYWJlbCA9ICJQYXVzZSIsIG1ldGhvZCA9ICJhbmltYXRlIiwKICAgICAgICAgICAgICAgYXJncyA9IGxpc3QoTlVMTCwgbGlzdChtb2RlID0gImltbWVkaWF0ZSIsIGZyYW1lID0gbGlzdChkdXJhdGlvbiA9IDApLCByZWRyYXcgPSBGQUxTRSkpKQogICAgICAgICkKICAgICAgKQogICAgKSwKICAgIHRpdGxlID0gIlRpbWU6IDAiCiAgKQoKIyBDb252ZXJ0IHRvIHBsb3RseSBvYmplY3Qgd2l0aCBmcmFtZSBpbmZvCnBiIDwtIHBsb3RseV9idWlsZChwKQoKIyBJbmplY3QgY3VzdG9tIHRpdGxlcyBpbnRvIGVhY2ggZnJhbWUKZm9yIChpIGluIHNlcV9hbG9uZyhwYiR4JGZyYW1lcykpIHsKICB0IDwtIHRpbWVfc2VxW2ldCiAgZXJyIDwtIHNpZ25pZihtYXhfZXJyb3JfYmV0d2Vlbl9ib3RoX2FwcHJveFtpXSwgNCkKICBwYiR4JGZyYW1lc1tbaV1dJGxheW91dCA8LSBsaXN0KHRpdGxlID0gcGFzdGUwKCJUaW1lOiAiLCB0LCAiIHwgRXJyb3I6ICIsIGVycikpCn0KYGBgCgpgYGB7cn0KZmlnICAjIERpc3BsYXkgdGhlIHBsb3QKYGBgCgoKYGBge3IsIGZpZy5oZWlnaHQgPSA2LCBvdXQud2lkdGggPSAiMTAwJSIsIGZpZy5jYXAgPSBjYXB0aW9uZXIoIkNhcHRpb24iKX0KcGIKYGBg